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We report on the observation of the Josephson effect between two strongly interacting fermionic 
superfluids coupled through a thin tunneling barrier. We prove that the relative population and 
phase are canonically conjugate dynamical variables, coherently oscillating throughout the entire 
crossover from molecular Bose-Einstein condensates (BEC) to Bardeen-Cooper-Schrieffer (BCS) 
superfluids. We measure the plasma frequency and we extract the Josephson coupling energy, both 
exhibiting a non-monotonic behavior with a maximum near the crossover regime. We also observe 
the transition from coherent to dissipative dynamics, which we directly ascribe to the propagation of 
vortices through the superfluid bulk. Our results highlight the robust nature of resonant superfluids, 
opening the door to the study of the dynamics of superfluid Fermi systems in the presence of strong 
correlations and fluctuations. 


The Josephson effect is a pristine example of a macro¬ 
scopic quantum phenomenon, disclosing the broken sym¬ 
metry associated with the superfluid state |Tj. On a very 
fundamental level, it allows to pinpoint the most elusive 
part of the superfluid order parameter, the phase, through 
a measurable quantity, a particle current [2]. Furthermore, 
being based on tunneling processes, Josephson dynamics 
provides fundamental insights into the microscopic proper¬ 
ties of superfluids and their robustness against dissipative 
phenomena [3]. Since its discovery, Josephson effect has 
been demonstrated for a variety of fermionic and bosonic 
systems EHH]. However, it has so far eluded observation 
in BEC-BCS crossover superfluids mm realized by ul¬ 
tracold Fermi gas mixtures close to a Feshbach resonance 
nans]. The interest in these systems is twofold: on the 
one hand, they encompass the two paradigmatic aspects 
of superfluidity within a single framework: Bose-Einstein 
condensation of tightly bound molecules and BCS super¬ 
fluidity of long-range fermion pairs |13j . Moreover, in the 
resonant regime where the pair size matches the interpar¬ 
ticle spacing, they exhibit universal properties, sharing 
analogies with other exotic strongly-correlated fermionic 
superfluids, from cuprate superconductors to nuclear and 
quark matter HZ] EH]. 

In this work, we report on the observation of the Joseph¬ 
son effect in ultracold gases of 6 Li atom pairs across the 
BEC-BCS crossover. Our Josephson junction consists of 
two superfluid reservoirs, weakly coupled through a thin 
tunneling barrier. For all interaction regimes, we detect 
coherent oscillations of both the pair population imbalance 
AN = Nl — Nr and the relative phase p = tl~Tr across 
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the junction, measured in situ and after time-of-flight ex¬ 
pansion respectively. We prove these two observables to be 
dynamically conjugate [2], directly unveiling macroscopic 
phase coherence in these strongly-correlated fermionic 
superfluids. We measure the plasma frequency wj of the 
oscillations around the crossover region, from which we 
extract the Josephson coupling energy Ej [3]. Both quan- 



FIG. 1. A: Sketch of the experimental apparatus. A Josephson 
junction is realized by bisecting trapped superfluids of 6 Li atom 
pairs with an optical barrier only a few times wider than the 
correlation length of the system. The dynamics is monitored 
by recording number imbalance and relative phase between the 
two reservoirs via in-situ (B) and time-of-flight (C) imaging, 
respectively. 

tities exhibit a non-monotonic behavior crossing over from 
the BEC to the BCS regime, reaching a maximum around 
the resonance region of unitary-limited interactions. For 
critical parameters of the barrier height or number im¬ 
balance, the coherent dynamics turns into a dissipative 
one. Throughout the whole BEC-BCS crossover, we find 
that the dissipative flow is triggered by the propagation 
of vortices in the superfluid bulk, which nucleate within 
the barrier region via phase slippage processes ehi [la¬ 
in our experiment, we produce fermionic superfluids 
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FIG. 2. Josephson oscillations of ultracold Fermi superfluids. A: Time evolution of z(t) (blue) and of ip(t) (red) in the BEC 
limit at 1/kFO, = 4.3 and Vo/y = 1.0. Each data point is the average of at least five independent measurements. Error bars are 
the corresponding standard deviation. Solid lines are fit to the data with a sinusoidal function. The oscillation frequencies of 
z(t) and of ip(t) obtained from the fit are 13.9(1) Hz and 13.8(2) Hz, respectively. The phase shift between the two oscillations 
is 1.1(1) 7r/2. B: Same as in A but at unitarity for Vo//u. = 1.1. Here the phase is detected after a rapid sweep to the BEC limit 
(see text). The extracted oscillation frequencies are 12.8(l)Hz and 12.6(3) Hz, respectively. The relative phase shift is equal to 
1.2(2) 7r/2. The shaded regions reflect the fit uncertainties. 


of about N = 10 5 atom pairs EDE2], confined in a har¬ 
monic potential of frequencies uj x = u>q = 2n x 15 Hz, 
uj y = 2n x 150 Hz and uj z =2ttx 170 Hz. The inter-atonric 
interaction, parameterized by the s-wave scattering length 
a , is finely tuned by exploiting a broad Feshbach reso¬ 
nance at 832 G [151115] ■ Since strong interactions easily 
foster the development of dissipative processes OHM], 
we engineer a thin optical barrier whose width is only 
a few times larger than the superfluid coherence length 
of the system [9]. This is realized by focusing onto the 
atomic cloud a strongly anisotropic laser beam at 532 
nm, blue-detuned with respect to the main optical transi¬ 
tion of Lithium atoms. At the trap center, the beam is 
Gaussian-shaped with 1/e 2 beam waist of 2.0(2) fin i and 
840(30) /mr along the x and y direction, respectively (see 
Fig. 04). The repulsive barrier is homogeneous along the 
y and z axes while bisecting the cloud on the weak i-axis 
in two reservoirs. We parametrize the barrier height with 
the potential peak Vo, felt by one atom pair, which is 
adjusted by controlling the power of the barrier beam. 
The dynamics of the system is then induced as it fol¬ 
lows: initially, the trap center is axially displaced with 
respect to the barrier position |22j . creating an initial 
non-zero population imbalance zq = AN/N between the 
two reservoirs. Here N = Nr + Nr is the sum of the pair 
populations of the two wells. A non-adiabatic movement 
of the trap center back onto the barrier position creates 
a non zero chemical potential difference Sfi o, triggering 
the superfluid dynamics. Initially, we focus our study 
on the regime of small excitations, working at the lowest 


detectable initial imbalance zq = 0.03(1) and for barrier 
heights Vo exceeding the bulk chemical potential fi, i. e. 
in the tunneling regime. Here the system Hamiltonian 
can be written in terms of only two macroscopic dynam¬ 
ically conjugate variables: the relative phase ip and the 
relative population AN/2 [2|- Quite generally, for small 
oscillations, the Hamiltonian reduces to the sum of two 
energy terms: Ejip 2 /2 and Ec/2(AN /2) 2 . The first one 
depends on the Josephson tunneling energy Ej and it 
favors the coherent flow of particles through the junction. 
The second term represents the charging energy, being 
Ec the energy cost to add a single pair in one reservoir 
HIS!. Consequently both imbalance and phase undergo 
harmonic oscillations, out-of-phase by 7t/2, at a plasma 
frequency that is independent from zq and given by: 

Hujj = \J EcEj ( 1 ) 

We confirm this expectation by studying the evolution 
of z(t) and <p{t) from absorption images, recorded in-situ 
and after time-of-flight expansion respectively, see Fig. 0 
To determine <p(t) also in the strongly interacting regime, 
we perform a 200 fis fast ramp to the BEC side of the 
resonance [32, 2S] ■ This reduces the detrimental effects 
of collisions during the expansion, that would completely 
wash out the visibility of the interferogram [25] , In Fig. 
0we present an example of the evolution of z(t ) and ip(t) 
for a molecular BEC and a unitary Fermi gas, respectively. 
In both cases, the two quantities oscillate at the same 
frequency, with a relative phase shift of 7r/2, within error 
bars. Our measurement proves the conjugate nature 
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of phase and number imbalance, representing a direct 
proof of the macroscopic phase coherence of these systems. 
Since this trend is confirmed also for a crossover superfluid, 
whose phase can be inferred only after ramping to the 
BEC limit, it indicates that the magnetic sweep does not 
appreciably distort the measurement of 

In order to understand how the Josephson dynam¬ 
ics is influenced by interactions through the BEC-BCS 
crossover, we extract the Josephson frequency loj by mea¬ 
suring z(t) for a fixed barrier height Vo = 1.2(1 )E F , 
at different values of the interaction strength that we 
parametrize via the dimensionless parameter l/k F a. 
Here k F is the Fermi wave-vector associated to a non- 




— 1 /k F a 

FIG. 3. Plasma frequency and Josephson coupling energy 
through the BEC-BCS crossover. A: ujj (blue circles) nor¬ 
malized to the bare trap frequency wo for Vo/Ef = 1-2(1) 
as a function of 1 /kFa. Each data point is the average of 
at least five independent measurements. Error bars are the 
corresponding standard deviation. B: Calculated charging 
energy Ec multiplied for the number of pairs N (red solid 
diamonds) and Josephson coupling energy Ej (blue circles) 
normalized to the Fermi energy E p extracted from Eq. 0. 
In both panels blue solid lines are a guide to the eye. 

interacting gas of N trapped fermions, whose Fermi energy 
is given by E F = = huj(6N) 1 / 3 (m is the mass of 

one 6 Li atom and w = (ui x uj v u) z ) 1 / 3 ). In Fig. [|A we 
present the behavior of wj normalized to u>o as a function 
of l/k F a. Notably, the extracted wj exhibits a non¬ 
monotonic evolution across the BEC-BCS crossover, with 


a maximum around the unitary limit. To gain further 
insights into this trend, by exploiting Eq. 0 we derive 
Ej from the measured values of the plasma frequency 
wj combined with the computed charging energy. Ec 
is related to the inverse of the system compressibility, 
and it has been derived from an extended Thomas-Fermi 
model, based on a generalized Gross-Pitaevskii equation 
for the pairs wavefunction that accounts for the correct 
chemical potential obtained by Monte Carlo calculations 
across the entire crossover |22j . In Fig. [3j3 we show Ej 
and N x Ec normalized to the Fermi energy E F as a 
function of l/k F a. While the charging energy increases 
monotonically moving from the BEC to the BCS regime, 
the Josephson coupling energy reflects the behavior of 
wj, reaching a maximum close to unitarity. Interestingly, 
our observation qualitatively reproduces the trend for 
the maximum Josephson current derived from numerical 
approaches [271 Efij which explicitly account for the com¬ 
posite nature of the superfluid pairs. For Vo /i, this 
behavior can be ascribed [28] to the competition of two 
different critical velocities at the crossover region, set by 
sound and pair-breaking excitations that are predominant 
on the BEC and BCS side respectively f29]. In turn, our 
observation in the tunneling regime Vq > fj, can be quali¬ 
tatively understood by noting that in both the BEC and 
BCS limits one expects 

Ej ~ JCN 0 . (2) 

Nq is the total number of condensed pairs, while K, is the 
tunneling term that depends on the chemical potential 
and the barrier properties [2J1 20]. In the deep BEC 
limit, almost all pairs are condensed, No — N and Ej ~ 
K.N. In the BCS limit Nq/N oc A/E f , where A is the 
superfluid gap, Eq. 0 reproducing the Ambegaokar- 
Baratoff formula at T = 0 [3]. Moving from the BEC to 
the unitary limit the increase of Ej reported in Fig. [3j3 
is associated with the growth of 1C. Indeed, in this regime 
N 0 decreases only slightly eh, while the increase of the 
interactions makes /r, hence /C, progressively larger. In 
contrast, towards the BCS limit, fi does not vary much 
with k F a, while N 0 is strongly decreased, causing a net 
reduction of Ej. Even if a microscopic derivation of 
the dependence of Ej upon Nq is presently missing in 
the strongly-interacting regime, the previous arguments 
suggest that a maximum of Ej should occur close to the 
unitary limit, in agreement with the experiment. 

So far, we have focused our studies on the regime of 
small excitations, characterizing the Josephson oscillations 
throughout the BEC-BCS crossover. Now, we investigate 
the system evolution once the dynamics is triggered by 
larger initial excitations. In our setup, this can be eas¬ 
ily done either by increasing the initial imbalance zq at 
fixed Vo, or vice versa. In this case, we observe a com¬ 
pletely different behavior, where the transport through 
the junction turns from coherent to dissipative. As an 
example, in Fig. [4] A we present the evolution of z(t) 
for z 0 = 0.04(1) and V 0 /E F = 2, at l/k F a = 0. This 
dynamics is in striking contrast with the one previously 
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FIG. 4. Onset of dissipative dynamics. A: Time evolution of z(t) at resonance for Vq/Ef = 2 and zo = 0.04(1) (blue diamonds). 
Solid line is an exponential fit to the data. Inset: Time evolution of ip(t) for the same initial parameters (green circles). The 
green solid line is a fit to a sawtooth function. Error bars are one standard deviation of five independent measurements. B: 
Vortex occurrence probability evaluated over a collection of 40 independent time-of-flight images (red diamonds, right axis), and 
normalized plasma frequency oij/wo (blue circles, left axis), as a function of Vq/Ef at l/fc_Fa = 0, for zo = 0.04(1). Error bars 
for the vortex occurrence are obtained by the Wilson score interval. Error bars for ojj/ljo represent the fit uncertainties. Inset: 
Typical image of one vortex, recorded after 10 ms time-of-flight expansion. C: Contour plot of u!j/qjq versus 1/k.Fa, for different 
Vq/Ef values indicated by the color code on the left. D: Same as in C for the probability of vortex occurrence. To speed up the 
data acquisition, zq = 0.12(2) and the statistics was performed over 20 images. In panels C and D the red line highlights the 
barrier value Vq/Ef = 1.2. 


reported in Fig. [2j the system is now characterized by 
a dissipative flow, which irreversibly tends to equilibrate 
the two reservoirs, similarly to what observed in Ref. [53] 
for crossover superfluids coupled through a mesoscopic 
channel. In this regime, it is interesting to study also 
the behavior of ip(t ). A priori, phase coherence between 
the two superfluids could be completely scrambled by the 
presence of dissipative processes. In turn, the visibility 
of the observed interference pattern remains high, and 
we are able to trace the initial evolution of see inset 
of Fig. [4] A. Also the dynamics of this observable signif¬ 
icantly differs from the one in Fig. [2j p(t) grows now 
linearly in time, at a rate much faster than ujq and com¬ 
parable with the initial chemical potential imbalance SfiQ- 
This behavior can be explained by noticing that, once the 
initial charging energy EqZqN 2 exceeds the Josephson 
coupling Ej , the phase increases as fup(t) ~ Ecz 0 Nt [53] , 
In the experiment, after an evolution time on the order of 
the axial trap period, despite the contrast of the interfero- 
grams remains constantly high in each image, we observe 


a significant increase of the shot-to-shot fluctuations that 
impedes to further follow a clear trend of ip(t). In order to 
understand how a slow resistive flow can coexist with the 
phase coherence between the two reservoirs, we investigate 
the microscopic origin of dissipation in our system. The 
high visibility of the interference pattern, detected at all 
evolution times, rules out pair-breaking processes, charac¬ 
terizing weakly interacting BCS superfluids, as the origin 
of the dissipative flow. Rather, our observation recalls the 
phenomenology typical of Helium systems, where resis¬ 
tive dynamics is established by phase-slippage processes 
and vortex nucleation. In fact, the entrance into this 
’’running phase” regime is expected to be accompanied 
by the nucleation of a vortex within the link region m- 
A complete phase slip ip = 2n occurs when the vortex 
annihilates within the barrier region. This gives rise to 
the macroscopic quantum self trapping of the relative 
population [221123] - However, depending on the system 
configuration, and for Vq ~ /i, such a topological defect 
may escape the low-density region mmm entering the 
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superfluid bulk, before annihilation. The propagation of 
the vortex through the superfluid bulk acts as a dissipa¬ 
tive channel that gives rise to a resistive flow that leads to 
an exponential decay of z(t). This mechanism can indeed 
occur in our crossover superfluids: the three-dimensional 
character of our junction, combined with the coupling to 
the transverse modes favored by the strong inter-particle 
interactions, may facilitate the leakage of vortices from 
the barrier region |32 . Indeed, by performing a statistical 
study over several time-of-flight images recorded after 
some time evolution in the trap }22] , we detect with non¬ 
zero probability the presence of topological defects, which 
appear as density depletions in the expanded clouds, see 
inset in Fig. [4] B. By measuring their oscillation period in 
the trap after switching off the barrier we identify them 
as solitonic vortices [22] [33]. The intimate connection be¬ 
tween the breakdown of the Josephson oscillations and the 
appearance of vortices is further confirmed by the study 
presented in Fig. [4]B. Here, we show the behavior of the 
Josephson frequency wj at unitarity as a function of Vo, 
together with the occurrence of defects collected for each 
Vo value over a statistical ensemble of 40 images. No vor¬ 
tices are detected until coherent oscillations are observed 
(Vo /Ep < 1.5), topological defects showing up only after 
the stop of coherent oscillations. Interestingly, the inter¬ 
connection between the quench of the coherent dynamics 
and the vortex nucleation is not peculiar of the unitary 
point, but it extends over the whole BEC-BCS crossover 
region. This can be clearly observed by comparing Fig. 
[4]C and D, where the measured uij is contrasted to the 
vortex occurrence probability, as a function of Vo/Ep and 
1/kpa. One can see how the trend of the first observable 
perfectly mirrors the behavior of the second one for all in¬ 
teraction regimes. In particular, Fig. [4] D highlights once 
more the particular robustness of the crossover superfluid, 


which shows a reluctance to the formation of topological 
defects while keeping the highest Josephson frequency. 
We also note that the results of our measurements differ 
from the ones reported when investigating the limit of 
vanishingly low barriers Vo <C pt, where plrononic excita¬ 
tions and pair-breaking effects, rather than vortices, cause 
the breakdown of superfluidity in the BEC and BCS side 
respectively [25] . 

In conclusion, we have observed coherent Josephson 
dynamics between two weakly-coupled superfluids 
of 1 2 3 4 5 6 7 8 * 10 Li atoms across the whole BEC-BCS crossover. 
The behavior of the Josephson coupling energy Ej, 
maximum nearby unitarity, together with the response 
to excitations, minimum at the resonance, provides 
further evidence of the high-Tq< nature of these fermionic 
superfluids also in the coherent tunneling regime. Our 
work paves the way to the study of the interplay between 
elementary and topological excitations in the dissipative 
dynamics obtained by varying the height and width 
of the interwell barrier and to the measurement of 
the superfluid gap in close analogy with tunneling 
experiments in superconductors m nsi- Moreover, 
extending our studies of the tunneling dynamics above 
the condensation temperature Tq may give insights 
on the role of phase fluctuations in the regime where 
preformed, non condensed pairs appear in the system m- 
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S.I. SAMPLE PREPARATION 

Fermi superfluids are prepared by evaporating a 
two-component mixture of the lowest hyperfine states 
of 6 Li in an optical dipole trap. We employ the 
|F = 1/2 ,mp = ±1/2) states (unless otherwise stated), 
labeled as |1) and |2). Following the procedure described 
in Ref. [ 33 ] , the atomic sample evaporated at a magnetic 
field of 832 G, on top of a broad Feshbach resonance, 
finally resulting in the production of a superfluid cloud 
of 1.0(1) xlO 5 atoms per spin state. At the end of the 
evaporation, the magnetic field is adiabatically ramped 
to the desired value, allowing fine tuning of the atomic 
scattering length a, evaluated using the magnetic field 
dependence a(B) reported in [ 55 ] . The trap is realized 
by two infrared laser beams, crossing with an angle of 
14°, in the a:y-plane. The primary beam has a wavelength 
of 1064 nm and it is focused on the atomic sample on a 
waist of 45 fi m. The secondary beam has a wavelength 
of 1070 nm and waists on the atoms of 100 and 45 /tin. 
The final harmonic potential is characterized by radial 
trapping frequencies of about 150 and 170 Hz along the y 
and z axes respectively and an axial frequency of about 
15 Hz. The magnetic field curvature of the Feshbach 
coils provides additional confining (anti-confining) po¬ 
tential along the x and y (z) directions, corresponding 
to trap (anti-trap) frequencies of uj x , y ~ 27 r x 8 Hz and 
w z = \' f 2ui X) y Hz at 832 G. These contributions sum in 
quadrature to the main optical trap frequencies, causing 
a variation of the overall trap frequencies of about 10% 
when spanning from the BEC to the BCS limit. 


A. Creation of the tunneling barrier 



FIG. SI. Equation of state of the unitary Fermi Gas. Main: 
Normalized compressibility (blue circles) k/ko versus normal¬ 
ized pressure P/Pq. k and P are scaled by the respective 
quantities Ko and Po at the same local density for a non¬ 
interacting Fermi gas. Inset: Specific heat per particle (red 
circles) CV/IV/cb as a function of the reduced temperature 
T/T f . 


beam. This reduces the effect of atom diffusion during 
the imaging pulse, which can distort the profile of the bar¬ 
rier beam carved in the atomic cloud. By extrapolating 
the barrier size to zero intensity we find uq=2.0(2) /rm. 
Numerical simulations of the dynamics in the BEC limit 
(see section S.V), in which the beam waist enters as a 


fixed parameter, further confirm the measured value of 


We produce the tunneling barrier by focusing on the 
atoms an anisotropic laser beam at 532 nm, by using a sin¬ 
gle aspheric lens (NA^ 0.6). This configuration produces 
a repulsive sheet of light which bisects the trapped cloud. 
The repulsive barrier is first characterized on a replica of 
the optical system implemented in the experiment. The 
measured beam waists, detected with a CCD camera, in 
the focus position are 1-9(1) g m and 840(30) /im, along 
the x and y- axis, respectively. To extrapolate the barrier 
width in the final set-up, we study the atomic density 
profiles in the BEC regime, exploiting the Thomas-Fermi 
approximation 77,(r)=max((/i— V(r))/g,0), where y is the 
chemical potential, V (r) is the double well-potential and 
g the interaction coefficient. In-situ absorption images of 
the atomic density distribution in the ajy-plane are taken 
by progressively decreasing the intensity level of the probe 


B. Control of the initial population imbalance and 
of the initial velocity 

The position of the 1070 nm secondary beam, rela¬ 
tive to the main one, is finely adjusted by changing the 
radio-frequency of the acousto-optic modulator control¬ 
ling its optical power. The change in the RF signal results 
in a shift of position of the beam focus. The 1070 nm 
beam mainly provides the trap axial confinement and its 
movement results only in an axial shift of the trap mini¬ 
mum. The barrier potential is then held fixed at a certain 
place which we identify as the final position of the trap 
minimum. Evaporation is performed in a trap displaced 
with respect to the barrier position in order to create 
a population imbalance among the left and right sides. 
The dynamics is triggered by non-adiabatically shifting 








the trap minimum over the barrier position. It is worth 
noticing that for the typical initial imbalances employed 
in this work, zq = 0.03(1), the maximum atomic velocity 
relative to the barrier is of the order of v = 0.03(1)uf, 
where up is the Fermi velocity. For any interaction regime 
this value lies well below the recently measured Landau 
critical velocity }3E[ for breaking the superfluid. 


S.II. THERMOMETRY AND EQUATION OF 
STATE 


A first indication about the temperature at the end 
of evaporation is obtained by measuring the condensed 
fraction in the BEG regime. Without the barrier, we esti¬ 
mate a condensed fraction always above 90% of the total 
pair number. The condensed fraction is determined by 
time of flight images recored after a Feshbach magnetic 
field sweep towards the BEC limit (690 G), by fitting 
the momentum distribution with a 2D bimodal function. 
Evaporation with the barrier on at any target height does 
not result in a change of the condensed fraction. More¬ 
over, the condensed fraction does not vary appreciably 
over the experimentally accessed evolution times, for all 
the different tunneling regimes. We further confirm the 
presence of a superfluid state at unitarity by determining 
the equation of state of the gas, following the procedure 
developed in m- In particular, we observe a peak in 
both the compressibility as a function of the pressure and 
in the specific heat, as a function of T/Tp, at a critical 
temperature T c = 0.15(3)Tf, where Tp = Ep/kp ((being 
kp the the Boltzmann constant), as shown in Fig. SI and 


in the inset, respectively. Our measured value of T c is 
in fair agreement with the one reported in m- We can 
extract the temperature of the gas by employing the virial 
expansion on the data. In this way, we estimate a degener¬ 
acy parameter at the trap center of T/Tp = 0.07(2). This 
gives a ratio T/T c below 0.6, which ensures a superfluid 
fraction around unity for the resonantly interacting Fermi 
gas [38] , 


S.III. PHASE MEASUREMENT 

We study the evolution of the relative phase by detect¬ 
ing the interference fringes arising from the two expanding 
clouds after a time of flight of typically 15 ms. The phase 
is extracted by fitting the resulting interferogram by a 2D 
Gaussian, modulated by a cosine function, of the form 
n(x, y) = Ae~ x ^ w ^e~ v / w y x (1 + Bcos(kx + ip)). In the 
BEC side of the resonance, we can extract the relative 
phase at the same magnetic field at which the dynamics 
of the imbalance is observed. For the acquisition of the 
data in the unitary limit and in the BCS side, we em¬ 
ploy a balanced mixture of the |1) — |3) states, where |3) 
corresponds to the |.F = 3/2, tof = —3/2) state at low 
magnetic fields. In this case, the sample is produced 
by performing the evaporation at 690 G, on top of the 



Time (ms) 

FIG. S2. Superfluid versus normal tunneling. A: The super¬ 
fluid oscillations (blue circles) are shown for different barrier 
heights. The values of Vo/e are: 0.25(2) (top), 1.50(3) (center), 
1.75(4) (bottom). B: The normal gas oscillations (red circles) 
are shown for different barrier heights. The values of Vo/e 
are 0.40(2) (top), 0.60(2) (center), 0.81(3) (bottom). In all 
the plots, the error bars are one standard deviation of five 
independent measurements. 


|1) — |3) Feshbach resonance. This resonance has a width 
significantly smaller than the one at 832 G. This facilitates 
the fast magnetic field sweep from the unitary to the BEC 
regime for increasing interference patterns visibility [39| . 
We have also checked the 11) — 13) population imbalance to 
evolve with the expected coupled dynamics with respect 
to the phase, and in accordance with the measurements 
taken on the 11) — 12) mixture in the same initial condition. 
We can also observe the transition towards completely 
disconnected superfluids when Vo > 4/i.The entrance in 
this disconnected regime is signaled by an increase in the 
variance A p 2 that approaches to 2-7T. 


S.IV. SUPERFLUID VERSUS NORMAL 
TUNNELING DYNAMICS 


In this section we compare the dynamics of a superfluid 
molecular BEC (1 /kpa = 4.6) with the one of a spin 
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polarized Fermi gas at T/Tp ~ 0.1. The latter is produced 
by evaporating a balanced |1) — |3) mixture at 300 G 
m- The scattering length is brought to zero by setting 
the Feshbach field around 570 G, in the vicinity of the 
zero-crossing of the |1) — |3) Feshbach resonance. The 
experimental results are shown in Fig. |S2| and Fig. |S3| 
Here, it is convenient to rescale the barrier height to the 
mean energy per particle e instead than the chemical 
potential /. i . In the BEC limit e = I/z, while for the 
normal Fermi gas e = | Ep. As already discussed in the 
main text, the superfluid undergoes coherent undamped 
oscillations also for barrier heights larger than e (see 
Fig. S2 A). For the normal, non interacting Fermi gas 
this is not the case. By increasing the barrier height, the 
dynamics of the system is characterized by an increasing 
damping of the oscillations, see Fig. S2 B). For Vo/e > 1, 
no oscillations are detected, the motion of the normal gas 
becoming overdamped in this regime. Moreover, while the 
superfluid systems exhibits a strong renormalization of 
the frequency as Vo is increased, the normal gas oscillates 
always at the bare trapping frequency wo, as shown in 

Fig. 1531 
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FIG. S3. Superfluid and normal gas oscillations versus the 
barrier height. Oscillation frequency, in units of the trap 
frequency u>o, as a function of Vo/e for the superfluid (blue cir¬ 
cles) and the normal gas (red diamonds). In all the plots, the 
vertical error bars are one standard deviation of five indepen¬ 
dent measurements, while the horizontal error bars take into 
account the experimental uncertainties of the atoms number 
and of the barrier width. 


S.V. EXTENDED THOMAS-FERMI MODEL 

In this section we briefly discuss the theoretical frame¬ 
work we use both for the determination of Ec in Fig. 3B 
of the main text, and for the numerical simulations of the 
system dynamics. In particular we make use of the so 
called extended Thomas-Fermi model (ETFM) pTTiTFfTTj . 
In the ETFM approach, the wavefunction >3/ obeys the 


equation 


ih 


d'F 

dt 


h 2 

= -V 2 4- + 2I/^ + 2F(2|^| 2 )T, (S.l) 

4m 


where |T(r, t)| 2 is the pair density normalized to the 
number of pairs N (the total number of particles being 
2 N), V(r) is the external potential felt by the atoms 
(having mass m). The nonlinear term in Eq. (S.l I is 
defined as F(n) = J^-, where n = 2| SR | 2 is the total 
fermion density and £ is the energy per particle, which 
has been determined across the BEC-BCS crossover by 
Monte Carlo calculations [44,. The choice of the nonlinear 
term F guarantees that the chemical potential of the bulk 
system, /i, is correct. By writing 2F(2|'F| 2 ) = /(|4/| 2 ) 
one has in the BEC side /(|SR| 2 ) = 47r |4>| 2 (where 

M = 2m is the mass of the pair and <zm= 0.6a is the 
molecular scattering length) and at unitarity /(|'F| 2 ) = 
2^(6 tt 2 ) 2 / 3 |4'| 4 / 3 (where £ is the Bertsch parameter: 
£ « 0.370). The potential 2H = V p acting on the pairs 
has the form 


V p (r) = \M (uj 2 x x 2 + oo 2 y y 2 + u 2 z z 2 ) + V 0 e ~ 2 *> 


(S.2) 


where u>i (i = x , y , z) are the trap frequencies and Vo an d 
w are the barrier height and width, respectively. In all 
numerical simulations the values of w, Vo, Wj, N and a 
are fixed to the experimental ones. The initial population 
imbalance is created, in the framework of the ETFM, as 
it follows: the center along x parabolic term in Eq. |S.2| 
is suddenly displaced to an initial value xq and then put 
back in x = 0. We considered different realistic times for 
this process and verified that the effect on the tunneling 
dynamics is negligible. From the solution 4/ (r, t) one can 
extract the number of pairs Np{t), Np(t) in the left and 
right well, respectively, and consequently z{t). 

The ETFM plasma frequency is obtained from the nu¬ 
merical solution of Eq. (S.l I with an initial imbalance 
z(t = 0) = zq —> 0: our results are plotted in Fig. S4 
of this supporting material as red squares. We point 
out that, while the theoretical results are obtained in 
the limit zq —> 0, the experimental values for wj are 
taken extracting the frequency of population oscillation 
using the lowest detectable initial imbalance z 0 = 0.03(1). 
However, the values of the plasma frequency predicted 
by the time-dependent ETFM are inside the error bars 
of the experimental data, from the BEC to the unitar¬ 
ity limit (see Fig. S4). The agreement is eventually 


lost in the BCS side where the pair breaking effects 
on the tunneling are not correctly taken into account 
by the ETFM. The charging energy Ec can be calcu¬ 
lated using the ETFM and it is given by Eq = 2U/N 
where U = N , the derivative being computed at 
Np = N/2. yi oc is the local chemical potential given 
by yioc = f dr ^ (V^> a ) 2 + V v </> 2 a + f{N a c/) 2 a )4> 2 a with 
a = L,R and (f> a are the left and right Wannier wave- 
functions computed from the ground and first excited 
states of the time-independent version of Eq. (S.l) [45] , 
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FIG. S4. Josephson frequency ojj in units of the trap frequency 
ujq for a fixed barrier height Vo/Ep = 1.2(1) as a function of 
1/fcFa around the resonance region: blue circles are the experi¬ 
mental findings, while red diamonds are the frequencies of the 
small oscillations calculated by numerically solving the time- 
dependent ETFM. The error bars in the experimental data 
are one standard deviation of five independent measurements. 


We studied the dependence of the chemical potential /x 
upon the barrier height verifying that the effect is of few 
percent (as expected since /x is a bulk property) and that 
an excellent estimate of IA in our system is obtained by 
directly using /x, computed with the barrier turned on, 
instead of using /x; oc . The result of our study is plotted 
in Fig. 3B in the main text as a red diamonds. Typical 
values of Vo for which the relative phase is well defined are 
Vo > /x, with /x ~ 100 hjjj x on the BEC side (1 /kpa ~ 4.5) 
and /.x ~ 500 hco x at the unitary limit. We investigated the 
values of Vo for which the two-mode approximation |46j 
gives a good description of the small oscillation frequency. 
In the two-mode model applied to Eq. (S.l), with the 
normalization of the wavefunction fixed to N to have the 
correct chemical potential, the tunneling parameter K. is 
the energy difference per particle between the first excited 
state and the ground state of the time-independent ETFM. 
By comparing the numerical results obtained from the 
numerical solution of Eq. (S.lI with the two-mode model 
findings we found that a good agreement is present for 
Vo > 1.4/x in the BEC limit and Vo > 1.2/x at unitarity. 
In our system, we have Vo ~ 1.1 /x in the BEC side and 
we expect to have additional contributions beyond the 
two mode approximation. We also explored the system 
dynamics, using the ETFM, for larger starting imbalances 
and compared the simulation results with the experiment. 
In Fig. |S5| we report the measured oscillation frequency 
for a molecular BEC (1 /kpa = 4.6), a unitary Fermi 
gas and a BCS superfluid (1 /kpa = —0.5) as a function 
of the barrier height. In all cases the initial imbalance 
is z 0 = 0.03(1). The theoretical oscillations frequencies, 
predicted by the ETFM, for the molecular BEC and the 
Unitary Fermi gas are also shown. In these simulations 


FIG. S5. Oscillation frequency uij in units of the trap frequency 
(Jo as a function of the normalized barrier height Vo/Ef, for 
a molecular BEC (red diamonds: 1/kFd = 4.6), a unitary 
Fermi gas (blue circles: 1/kpa = 0) and a BCS superfluid 
(green squares:l/fcFa = —0.5). The initial imbalance is zo = 
0.03(1). The corresponding values of ETFM predictions for 
the BEC (red solid curve) and unitary gas (blue solid curve) 
are also shown. In the numerical simulations the starting 
imbalance is set to zo = 0.03 in agreement with the experiment. 
The vertical error bars are one standard deviation of five 
independent measurements, while the horizontal error bars 
take into account the experimental uncertainties of the atoms 
number and of the barrier width. 


we fix the value of zq to the experimental one. We observe 
that, in all the explored interaction regimes, increasing Vo 
leads to a progressive lowering of the oscillation frequency. 
As Vo exceeds a critical value V^ cr \ the oscillations of 
z(t) are eventually lost. In the BEC side one sees that the 
ETFM gives an excellent description also of the large os¬ 
cillation regime. This is not the case for the unitary Fermi 
gas, where the ETFM severely underestimates Vq CT> . even 
though gives a good description of the small oscillations 
as seen in Fig. [S4j 


S.VI. OBSERVATION AND DYNAMICS OF 
TOPOLOGICAL DEFECTS 


At a critical value of the initial population imbalance 
and/or height of the barrier, the system enters in a “run¬ 
ning phase” regime. This is accompanied by a continuous 
creation and annihilation of topological excitations inside 
the barrier region m Depending on the geometry of 
the system and on the strength of the inter-particle in¬ 
teraction 0UE9) the topological excitations can leave 
the low density region and escape inside the bulk of the 
system where they can be observed experimentally. The 
emergence of such dissipation mechanism is confirmed by 
the observation of vortices created at the center of the 
barrier and propagating in the direction of the superfluid 
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FIG. S6. Experimental procedure for detecting the topological 
defects. 


flow. Their presence is revealed by following the proce¬ 
dure schematically illustrated in Fig. [S6] We initialize 
the tunneling dynamics and we let the system evolve for 
50 ms, a time scale longer than half of the axial trapping 
period. Successively the barrier is slowly switched off by 
linearly decreasing the laser power from the initial value 
to zero in 80 ms. After a variable holding time, we turn 
off the optical confinement, and let the system expand for 
typically 10 ms after which an image of the atomic density 
distribution is recorded. Once entering the running phase 
regime, we typically observe with non-zero probability 
a single defect, appearing as a density depletion in the 
cloud, see Fig. 4A in the main text. Before the time of 
flight expansion, the scattering length is tuned to the 
BEC side of the Feshbach resonance (690 G) by applying 
a 10 ms linear ramp in the magnetic field. This procedure 
converts fermionic pairs into tightly bound molecules, 
emptying out the vortex core and enhancing the visibility 
of the topological defect j5T|. To gain more insights into 
the nature of the observed defects, we measure their axial 
oscillation period in the trap, once the barrier is removed. 
Fig. |S7| shows the vortex period Tjj, normalized to the 
trap period To=2n/ ujq , as a function of — 1/fc^a. In the 
inset we instead present the evolution of the position of 
the vortex in the trap for three different values of the 
interaction parameter. We note that the vortex period Tjj 


increases of about one order of magnitude moving from 
the BEC to the BCS side. This result is in agreement 
with the recent observation of solitonic vortices in BEC- 
BCS crossover superfluids, deterministically generated by 
phase-imprinting techniques ED- 

We can infer that, in our system, the phase slippage pro¬ 
cess is induced by a vortex ring formed at the low density 
barrier region in a plane perpendicular to the flow |48j . 
However, this kind of defect is unstable when it moves 
in an anisotropic trap ( ui v ^uj z ), as in our case, hence the 
ring will break into two vortex lines and finally will decay 
into a more stable structure which we can identify with a 
solitonic vortex in accordance with the large trap aspect 
ratio [52, 53. . In future work, it would be interesting to 
further investigate the mechanism underlying the dissi¬ 
pation of the superfluid flow due to the propagation of 



FIG. S7. Motion of topological defects. Main: Normal¬ 
ized period To/To of a single defect as a function of —X/kpa 
(blue squares), after removing the barrier. The blue line is a 
guide to the eye. The error bars correspond to the fit uncer¬ 
tainty. Inset: Axial position of a single defect versus time for 
l/fcFft=4.25, 1, 0 (red circles, blue triangles, and green squares, 
respectively). Here, the error bars are one standard deviation 
of five independent measurements in which the vortex was 
detected. 


topological defects as well as their dynamical evolution. 
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